clear
set more off

local inputs_dir "C:\Users\Administrator\Dropbox\Bond.Risk.Premia.Mitra.Xu\Data\combined_results_yu\data"
local outputs_dir "C:\Users\Administrator\Dropbox\Bond.Risk.Premia.Mitra.Xu\Data\combined_results_yu\outputs\outputs_replicate_fb"

* local SAMPLE_START_YEAR 1953
local SAMPLE_START_YEAR 1964
local SAMPLE_END_YEAR 2016
* local SAMPLE_END_YEAR 2007

local NEWEY_WEST_LAGS 12
local NEWEY_WEST_LAGS_QTR 4

local LAGN 1

/* ----------------------------------------------------------------------------------------------------------------------------------------------------------
1. set up data
---------------------------------------------------------------------------------------------------------------------------------------------------------- */

cd "`outputs_dir'"

* Fama-Bliss ZCB prices

use "`inputs_dir'\fama_bliss_prices.dta"
* format: Date_yyyymmdd (end of month price), price_1y, ..., price_5y
* note: prices are for a face value of 100

* convert date to datem
gen D = mod(Date_yyyymmdd,100)
gen M = mod((Date_yyyymmdd-D)/100,100)
gen Y = (Date_yyyymmdd - 100*M - D)/10000
gen datem = ym(Y,M)
format datem %tm

tsset datem

keep datem price_1y price_2y price_3y price_4y price_5y
order datem price_1y price_2y price_3y price_4y price_5y

/* ----------------------------------------------------------------------------------------------------------------------------------------------------------
2. Set up variables
---------------------------------------------------------------------------------------------------------------------------------------------------------- */

* log prices
foreach n of numlist 1/5{
	
	gen log_price_`n'y = log(price_`n'y/100)
	
}

* short rate
gen short = -100*log_price_1y

* ytm
foreach n of numlist 2/5{
	gen ytm_`n'y = -100*log_price_`n'y/`n'
}

* forward rates
foreach n of numlist 2/5{
	
	local lag_n = `n' - 1
	gen fwd_`n'y = 100*(log_price_`lag_n'y - log_price_`n'y)
	gen fwd_`n'y_spread = fwd_`n'y - short
	
}

* holding period return (annual)
foreach n of numlist 2/5{
	
	local lag_n = `n' - 1
	gen ret_`n'y = 100*(log_price_`lag_n'y - L12.log_price_`n'y)
	gen exret_`n'y = ret_`n'y - L12.short
	
}

gen ave_exret = (exret_2y + exret_3y + exret_4y + exret_5y)/4

gen slope = ytm_5y - short

save temp_dat.dta, replace

/* ----------------------------------------------------------------------------------------------------------------------------------------------------------
3. Monthly Return predictive regressions
---------------------------------------------------------------------------------------------------------------------------------------------------------- */

use temp_dat.dta, clear

local varlist1 ave_exret exret_2y exret_3y exret_4y exret_5y short fwd_2y fwd_3y fwd_4y fwd_5y slope

eststo clear
estpost summarize `varlist1' if datem>=ym(`SAMPLE_START_YEAR',1) & datem<=ym(`SAMPLE_END_YEAR',12), detail
esttab using "`outputs_dir'\replicate_fb_`SAMPLE_START_YEAR'_`SAMPLE_END_YEAR'.tex", replace cells("count mean sd min p25 p50 p75 max") nomtitle nonumber
esttab using "`outputs_dir'\replicate_fb_`SAMPLE_START_YEAR'_`SAMPLE_END_YEAR'.csv", replace cells("count mean sd min p25 p50 p75 max") nomtitle nonumber
eststo clear

* correlation table
eststo clear
estpost corr `varlist1' if datem>=ym(`SAMPLE_START_YEAR',1) & datem<=ym(`SAMPLE_END_YEAR',12), matrix
esttab using "`outputs_dir'\corr_table_`SAMPLE_START_YEAR'_`SAMPLE_END_YEAR'.tex", unstack b(2) p(2) noobs nostar note("Note: P-values in parentheses") title(Correlation matrix) nonumber nogaps replace compress
esttab using "`outputs_dir'\corr_table_`SAMPLE_START_YEAR'_`SAMPLE_END_YEAR'.csv", unstack b(2) p(2) noobs nostar note("Note: P-values in parentheses") title(Correlation matrix) nonumber nogaps replace compress substitute({table} {sidewaystable})
eststo clear

* predictive regressions
eststo clear
foreach n of numlist 2/5{
	qui regress exret_`n'y L12.fwd_`n'y_spread if datem>=ym(`SAMPLE_START_YEAR',1) & datem<=ym(`SAMPLE_END_YEAR',12)
	local r2_temp = round(`e(r2)',.001)
	eststo: qui newey exret_`n'y L12.fwd_`n'y_spread if datem>=ym(`SAMPLE_START_YEAR',1) & datem<=ym(`SAMPLE_END_YEAR',12), lag(`NEWEY_WEST_LAGS')
	qui estadd local sample_year "`SAMPLE_START_YEAR'-`SAMPLE_END_YEAR'"
	qui estadd local display_NWlags "`NEWEY_WEST_LAGS'"
	qui estadd local r_squared "`r2_temp'"
}
esttab using "`outputs_dir'\\fb_predict_ret_`SAMPLE_START_YEAR'_`SAMPLE_END_YEAR'.tex", replace scalars("sample_year Sample" "display_NWlags NW-lags" "r_squared R2") star(+ 0.10 * 0.05 ** 0.01 *** 0.001)
esttab using "`outputs_dir'\\fb_predict_ret_`SAMPLE_START_YEAR'_`SAMPLE_END_YEAR'.csv", replace scalars("sample_year Sample" "display_NWlags NW-lags" "r_squared R2") star(+ 0.10 * 0.05 ** 0.01 *** 0.001)
eststo clear

eststo clear
foreach n of numlist 2/5{
	qui regress exret_`n'y L12.slope if datem>=ym(`SAMPLE_START_YEAR',1) & datem<=ym(`SAMPLE_END_YEAR',12)
	local r2_temp = round(`e(r2)',.001)
	eststo: qui newey exret_`n'y L12.slope if datem>=ym(`SAMPLE_START_YEAR',1) & datem<=ym(`SAMPLE_END_YEAR',12), lag(`NEWEY_WEST_LAGS')
	qui estadd local sample_year "`SAMPLE_START_YEAR'-`SAMPLE_END_YEAR'"
	qui estadd local display_NWlags "`NEWEY_WEST_LAGS'"
	qui estadd local r_squared "`r2_temp'"
}
esttab using "`outputs_dir'\\slope_predict_ret_`SAMPLE_START_YEAR'_`SAMPLE_END_YEAR'.tex", replace scalars("sample_year Sample" "display_NWlags NW-lags" "r_squared R2") star(+ 0.10 * 0.05 ** 0.01 *** 0.001)
esttab using "`outputs_dir'\\slope_predict_ret_`SAMPLE_START_YEAR'_`SAMPLE_END_YEAR'.csv", replace scalars("sample_year Sample" "display_NWlags NW-lags" "r_squared R2") star(+ 0.10 * 0.05 ** 0.01 *** 0.001)
eststo clear
